%-------------------------------------------------------------------------
% Subroutine: Computation Historical Decomposition
%
% Estimation Codes for State-Level Results in  
% Ben-David, Itzhak, Sebastian Weber, and Pascal Towbin 
% "Inferring Expectations from Observables: Evidence from the Housing Market" 
% Review of Economics and Statistics
% ------------------------------------------------------------------------ 

%% Historical Decomposition
n_decomp=n_price_n;     % House price set as variable to decompose

for is=1:51
if is==37
else
ct=is;
state=State(t*ct,1);
state_vec(ct,1)=state;
[nstr, TT, ndraws,xx]=size(strerr);        
histirf_state(:,:,:,:)=zeros(n_var,TT,nstr,tel);
construct_y_sta;

J = zeros(n_var,n_var*p); J(:,1:n_var) = eye(n_var);

for tel3=1:tel
   AA=eye(size(A,1));
   start=1;
   constant=0;
 for tt=1:TT-start+1
     
     constant=constant+   J* AA*J'*gamma_g(n_var*p+1,:,ct)'; 
        AA=A(:,:,ct)*AA;

for strshoc=1:nstr;

    helpirf=squeeze(irfstr(1:n_var,strshoc,:,tel3,ct)) ;
    helpirf2=[zeros(n_var,tt-1),   helpirf(:,1:TT-tt+1).*repmat(strerr(strshoc,tt+start-1,tel3,ct),n_var,TT-tt+1) ];
    histirf_state(:,:,strshoc,tel3)= histirf_state(:,:,strshoc,tel3)+ helpirf2;
    histempt_state(:,tt,strshoc,tel3)=zeros(n_var,1)+constant; 
    
     for ll=1:p;
     histempt_state(:,tt,strshoc,tel3)= histempt_state(:,tt,strshoc,tel3)+AA(1:n_var,(ll-1)*n_var+1:ll*n_var)*y(start+p-ll,:)' ;
     end
     
   end
   end
end


histirf_temp=histirf_state;
y_hist=y(start+p:end,:);
y_hist(:,cum )=cumsum(y(start+p:end,cum),1);
histirf_state(cum,:,:,:)=cumsum(histirf_temp(cum,:,:,:),2);
histempt_state(cum,:,:,:)=cumsum(histempt_state(cum,:,:,:),2);

   
realdata = y_hist(:,n_decomp);



coll_histirf(:,:,:,:,ct)=histirf_state(:,:,:,:);
coll_hist_mean(:,:,:,ct)=mean(histirf_state(:,:,:,:),4);
coll_histdet_mean(:,:,:,ct)=mean(histempt_state(:,:,:,:),4);
coll_realdata(:,:,ct)=y_hist';
end
end



coll_histirf(:,:,:,:,37)=[];
coll_hist_mean(:,:,:,37)=[];
coll_histdet_mean(:,:,:,37)=[];
coll_realdata(:,:,37)=[];
state_vec(37,:)=[];
agg_realdata=mean(coll_realdata,3);
[max_hp]=squeeze(max(coll_realdata(n_decomp,:,:),[],2));


histempt_cent=mean(coll_hist_mean(:,:,:,:),4);
realdata=agg_realdata(n_decomp,:)';


histirf=coll_hist_mean;
histempt=coll_histdet_mean;

%%%%%%%%%%%%%%%%%%%%%%%%%%%Normalisierung
datelow=find(v_year==1996)-p;


norms=[n_price_n]; %normalize for prices only
 for ii=1:size(norms,2)
     for ihist=1:size(histirf,2)
     histirf2(norms(ii),ihist,:,:)=histirf(norms(ii),ihist,:,:)-histirf(norms(ii),datelow,:,:);
     histempt2(norms(ii),ihist,:,:)=histempt(norms(ii),ihist,:,:)-histempt(norms(ii),datelow,:,:);
     end
 coll_realdata(norms(ii),:,:)=squeeze(coll_realdata(norms(ii),:,:))-ones(size(coll_realdata(norms(ii),:,:),2),1)*squeeze(coll_realdata(norms(ii),datelow,:))';
 y_hist(:,norms(ii))=y_hist(:,norms(ii))-ones(size(y_hist,1),1)*y_hist(datelow,norms(ii));
 histirf(norms(ii),:,:,:)=histirf2(norms(ii),:,:,:);
 histempt(norms(ii),:,:,:)=histempt2(norms(ii),:,:,:);
 end
histirf_med=median(histirf(:,:,:,:),4);
histempt_med=median(histempt(:,:,:,:),4);